Guatemala Forest Fire Project

Project background

Research Questions:

1) When were fires most abundant in Guatemala?

2) Where are the firest most common in Guatemala?

3) What might be causing this surge in fires?

4) How do the fires affect biodiversity?

5) Is there an area at greater risk of fires and forest loss in Guatemala?

Data Acquisition and Organization

Fire data was acquired for 2019 from NASA. Specifically, data comes from the Near real-time (NRT) Suomi National Polar-orbiting Partnership (Suomi NPP) Visible Infrared Imaging Radiometer Suite (VIIRS) Active Fire detection product. For more information, please refer to NASA's Earthdata Website: https://earthdata.nasa.gov/earth-observation-data/near-real-time/firms/v1-vnp14imgt#ed-viirs-375m-attributes.

VIIRS (S-NPP) I Band 375 m Active Fire Product NRT (Vector data) data for 2019 was downloaded as a shapefile and visualized in ESRI ArcMap 10.8.1. Forest data comes from the National Institute of Forests in Guatemala (https://data.globalforestwatch.org/datasets/7935041390964af0a09763ff83c30b0e).

Using the VIIRS point data, I extracted values from the Forest Cover raster to determine the kind of forest affected by fires in 2019. Data from this shapefile was then exported as a CSV for ease of access.

# Import libraries
import pandas as pd
import plotly.express as px
import seaborn as sns
from matplotlib import pyplot as plt
import numpy as np
import geopandas as gpd 

# Load the Forest Fire data as a DataFrame, labelled VIIRS_fire_df
# I'm using the VIIRS (S-NPP) I Band 375 m Active Fire Product NRT (Vector data). 

VIIRS_fire_df = pd.read_csv(r"C:\Users\omara\OneDrive\Coding 2020 Resources\Guatemala_Forest_Project\Fires_forest_data.csv")
VIIRS_fire_df.head()

	FID	latitude	longitude	bright_ti4	scan	track	acq_date	acq_time	satellite	instrument	confidence	version	bright_ti5	FRP	daynight	fire_type	forest_type
0	275	13.978041	-91.286400	349.85	0.41	0.37	1/1/2019	1929	N	VIIRS	n	1	302.35	10.28	D	Vegetation Fire	No Forest
1	276	13.975160	-91.282066	340.52	0.41	0.37	1/1/2019	1929	N	VIIRS	n	1	308.64	49.20	D	Vegetation Fire	No Forest
2	277	13.974609	-91.285889	367.00	0.41	0.37	1/1/2019	1929	N	VIIRS	h	1	334.59	140.68	D	Vegetation Fire	No Forest
3	622	13.957825	-90.720909	342.94	0.43	0.38	1/1/2019	1929	N	VIIRS	n	1	304.29	5.31	D	Vegetation Fire	No Forest
4	659	13.977004	-90.686859	334.27	0.43	0.38	1/1/2019	1929	N	VIIRS	n	1	303.39	6.29	D	Vegetation Fire	Bosque Latifoliado Ralo

# Add Columns to the data frame to organize the data by month.

VIIRS_fire_df['Date'] = pd.to_datetime(VIIRS_fire_df['acq_date'])
VIIRS_fire_df['Month'] = VIIRS_fire_df['Date'].dt.month
VIIRS_fire_df['Month-str'] = VIIRS_fire_df['Date'].dt.strftime('%b')
VIIRS_fire_df['Month-str-full'] = VIIRS_fire_df['Date'].dt.strftime('%B')

VIIRS_fire_df.head()

	FID	latitude	longitude	bright_ti4	scan	track	acq_date	acq_time	satellite	instrument	...	version	bright_ti5	FRP	daynight	fire_type	forest_type	Date	Month	Month-str	Month-str-full
0	275	13.978041	-91.286400	349.85	0.41	0.37	1/1/2019	1929	N	VIIRS	...	1	302.35	10.28	D	Vegetation Fire	No Forest	2019-01-01	1	Jan	January
1	276	13.975160	-91.282066	340.52	0.41	0.37	1/1/2019	1929	N	VIIRS	...	1	308.64	49.20	D	Vegetation Fire	No Forest	2019-01-01	1	Jan	January
2	277	13.974609	-91.285889	367.00	0.41	0.37	1/1/2019	1929	N	VIIRS	...	1	334.59	140.68	D	Vegetation Fire	No Forest	2019-01-01	1	Jan	January
3	622	13.957825	-90.720909	342.94	0.43	0.38	1/1/2019	1929	N	VIIRS	...	1	304.29	5.31	D	Vegetation Fire	No Forest	2019-01-01	1	Jan	January
4	659	13.977004	-90.686859	334.27	0.43	0.38	1/1/2019	1929	N	VIIRS	...	1	303.39	6.29	D	Vegetation Fire	Bosque Latifoliado Ralo	2019-01-01	1	Jan	January

5 rows × 21 columns

# Subdividing data to get the number of fires by month. 
fire_data_by_month = VIIRS_fire_df.groupby(['Month', 'Month-str']).FID.count().reset_index()
fire_data_by_month = fire_data_by_month.rename(columns={'FID':'Count'})

fire_data_by_month

	Month	Month-str	Count
0	1	Jan	1688
1	2	Feb	2537
2	3	Mar	6023
3	4	Apr	14015
4	5	May	9144
5	6	Jun	1091
6	7	Jul	622
7	8	Aug	1010
8	9	Sep	471
9	10	Oct	461
10	11	Nov	816
11	12	Dec	1322

# Graph showing when fires were most abundant in Guatemala in 2019
fig1 = plt.figure(figsize=(10,8))
sns.set_style('darkgrid')
ax = sns.barplot(data=fire_data_by_month, x='Month-str', y='Count', palette=['#ffffb2','#fecc5c','#fd8d3c','#e31a1c'])
ax.set_title('Distribution of Fires in Guatemala by Month')
ax.set_ylabel('Number of Fires')
ax.set_xlabel('Months')
plt.show()

# Violin plot showcasing the distriubtion of VIIRS I-4 channel brightness temperature of the fire pixel measured in Kelvin.
# This gives us an idea of how hot the fires were burning when data was recorded at the sensor. 

fig2 = plt.figure(figsize=(10,8))
ax2 = sns.violinplot(data=VIIRS_fire_df, x='fire_type', y='bright_ti4')
ax2.set_title('Distribution of Fires in Guatemala by Type')
ax2.set_ylabel('Bright ti4')
ax2.set_xlabel('Fire Type')
plt.show()

On Violin Plots:

There are two KDE plots that are symmetrical along the center line. A white dot represents the median. The thick black line in the center of each violin represents the interquartile range. The lines that extend from the center are the confidence intervals. The violin plot also displays the 95% confidence interval.

#Visualize the Data using Mapbox Scattermap. I'm using a basic base map I made in Mapbox that highlights parks and forests. 
map = px.scatter_mapbox(VIIRS_fire_df, 
                        lat= "latitude", 
                        lon= "longitude", 
                        size='FRP', 
                        color= 'fire_type',
                        hover_name="bright_ti4", 
                        hover_data=['bright_ti4', 'fire_type'],
                        zoom=6.5,
                        size_max=50, 
                        title= 'Fires in Guatemala by Month 2019', 
                        labels={'fire_type': 'Fires'},
                       animation_frame='Month-str')


map.update_layout(mapbox_style="mapbox://styles/omarandres27/ckgbbi1dl441l19o6vbqubc6x", 
                mapbox_accesstoken='pk.eyJ1Ijoib21hcmFuZHJlczI3IiwiYSI6ImNrZjBhYm01czBkc3UyeHFpMjFtdXZ0cWUifQ.TpHxJmlbZoRC45LaYjmcSw')
map.update_layout(margin={"r":1,"t":30,"l":1,"b":0})

# Show our interactive map
# map.write_html(r"C:\Users\omara\OneDrive\Coding 2020 Resources\Guatemala_Forest_Project\Guatemala_Forest_Fire_Map.html", auto_open=True)
map.show()

from IPython.display import Image
img_bytes = map.to_image(format="png")
Image(img_bytes)

Now that we know when fires are most prevalent, let's figure out what kind of forest cover is most affected.

As mentioned above, forest data comes from the National Institute of Forests in Guatemala (https://data.globalforestwatch.org/datasets/7935041390964af0a09763ff83c30b0e). Researchers classified 308 high resolution RapidEye data to create 16 classes of forests for the year 2012.

Using ArcGIS ArcMap, I extracted the raster values representing different forest cover where they interestect with different points representing fires. This allows us to approximate the kind of forest that was predominantly affected by fires in 2019. From the data, we know that the Bosque Latifoliado Denso (Dense Broadleaf Forest) is the most affected by fires in 2019. This is also the type of forest cover most typical in the lowlands of northern Guatemala (for more information see: http://www.reddccadgiz.org/documentos/doc_1170376601.pdf).

# Subdividing data by forest type. What kind of Forest cover was most susceptible to fires in 2019. 

forest_data = VIIRS_fire_df.groupby('forest_type').FID.count().reset_index()
forest_data = forest_data.rename(columns={'FID':'Count'})
forest_data = forest_data.rename(columns={'forest_type':'Forest Cover'})
forest_data

	Forest Cover	Count
0	Bosque Conífero Denso	540
1	Bosque Conífero Ralo	742
2	Bosque Latifoliado Denso	17175
3	Bosque Latifoliado Ralo	7572
4	Bosque Mixto Denso	2005
5	Bosque Mixto Ralo	1373
6	Mangle Blanco Denso	50
7	Mangle Blanco Ralo	6
8	Mangle Negro Denso	9
9	Mangle Rojo Denso	142
10	No Forest	9586

# Graph of the forest cover most affected by fires in 2019.
fig4 = plt.figure(figsize=(10,8))
sns.set_style('darkgrid')
ax2 = sns.barplot(data=forest_data, x='Forest Cover', y='Count', palette=['#ffffcc','#c2e699','#78c679','#238443'])
ax2.set_title('Forest affected by fires')
ax2.set_ylabel('Amount of Fires')
ax2.set_xticklabels(ax2.get_xticklabels(), rotation=45,
    horizontalalignment='right')
plt.show()

Fires and Forests by Departments in Guatemala.

The following section will evaluate the amount of fires in 2019 by department, as well as Forest Cover Loss data from 2001 to 2019. From the data, is clear that the department of the Peten has been the center of fires in 2019. Unfortunately, the Peten has also seen the most Forest Canopy Loss since 2001.

This data was acquired using ArcGIs Pro's Aggregate Point Function. I aggregated the count of fires by department. This data was saved into te fire_aggregate_data_by_dept shapefile and used a geo data frame using Geopandas.

# Number of fires in 2019 by Department.  
geo_dept_fires = gpd.read_file(r"C:\Users\omara\OneDrive\Coding 2020 Resources\Guatemala_Forest_Project\Fire Data\fire_aggregate_data_by_dept.shp")
geo_dept_fires.head()

# Number of fires in 2019 by departments in Guatemala.
fires_by_dept = plt.figure(figsize=(10,8))
sns.set_style('darkgrid')
ax3 = sns.barplot(data=geo_dept_fires, x='ADM1_REF_1', y='COUNT_FRP', palette=['#9e0142','#d53e4f','#f46d43','#fdae61','#fee08b','#ffffbf',\
                                                             '#e6f598','#abdda4','#66c2a5','#3288bd','#5e4fa2'])
ax3.set_title('Number of Fires in 2019 by Departments in Guatemala')
ax3.set_ylabel('Number of Fires')
ax3.set_xlabel('Department')
ax3.set_xticklabels(ax3.get_xticklabels(), rotation=45,
    horizontalalignment='right')
plt.show()

# Choropleth map showing the regions in Guatemala that have been disproportionally affected by fires in 2019.
fire_dept_map = px.choropleth(geo_dept_fires, 
                          geojson=geo_dept_fires.geometry, 
                          locations=geo_dept_fires.index,
                          hover_name= 'ADM1_REF_1',
                          projection="mercator",
                          color='COUNT_FRP',
                          color_continuous_scale='Oryel',
                          labels={'COUNT_FRP': 'Number of Fires'},
                          title= 'Number of Fires in Guatemala by Department 2019')
fire_dept_map.update_geos(fitbounds="locations", visible=False)
fire_dept_map.update_layout(margin={"r":0,"t":36,"l":0,"b":0})
fire_dept_map.show()

# Show static version of the map
depat_fire_map = fire_dept_map.to_image(format="png")
Image(depat_fire_map)

How does the location of fires relate to Forest Cover Loss in the last 19 years?

Data for Forest Cover Loss comes from the University of Maryland Global Forest Change Project: http://earthenginepartners.appspot.com/google.com/science-2013-global-forest

‘Forest Cover Loss’ is defined as a stand-replacement disturbance, or a change from a forest to non-forest state, during the period 2000–2019.

Data was organized in the GTM_dept_ForestLoss shapefile and used as a geo data frame using Geopandas. The SUM column represents the sum of pixels (from 0 to 19) where there has been forest cover loss.

# Sum of Forest Canopy Loss betwween 2001 and 2019 by Departments in Guatemala. 

# Read the GTM_dept_ForestLoss shapefile that has the polygon geometry as well as the Forest Cover Loss data.
geo_df = gpd.read_file(r"C:\Users\omara\OneDrive\Coding 2020 Resources\Guatemala_Forest_Project\Forest Data\GTM_dept_ForestLoss.shp")
geo_df.head()

# graph showing distribution of areas with the most forest loss since 2001 (from the University of Maryland Project)

forest_loss = plt.figure(figsize=(10,8))
sns.set_style('darkgrid')
ax4 = sns.barplot(data=geo_df, x='ADM1_ES', y='SUM', palette=['#9e0142','#d53e4f','#f46d43','#fdae61','#fee08b','#ffffbf',\
                                                             '#e6f598','#abdda4','#66c2a5','#3288bd','#5e4fa2'])
ax4.set_title('Forest Cover Loss in Guatemala by Department 2001-2019')
ax4.set_ylabel('Forest Loss')
ax4.set_xlabel('Department')
ax4.set_xticklabels(ax4.get_xticklabels(), rotation=45,
    horizontalalignment='right')
plt.show()

# Choropleth map showing Forest Canopy Loss in Guatemala between 2001 and 2019. 
guate_map = px.choropleth(geo_df, 
                          geojson=geo_df.geometry, 
                          locations=geo_df.index,
                          hover_name= 'ADM1_REF_1',
                          projection="mercator",
                          color='SUM',
                          color_continuous_scale='Sunsetdark',
                          labels={'SUM': 'Forest Cover Loss'},
                          title= '2001-2019 Forest Cover Loss in Guatemala by Department')
guate_map.update_geos(fitbounds="locations", visible=False)
guate_map.update_layout(margin={"r":0,"t":36,"l":0,"b":0})
guate_map.show()

# Show static version of the map
depat_forest_map = guate_map.to_image(format="png")
Image(depat_forest_map)

Conclusions:

Forest Fire Data from SIFGUA: http://www.sifgua.org.gt/Incendio.aspx

print('Script Complete')

Script Complete